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Abstract. While it is generally thought that molecular 
outflows from young stellar objects (YSOs) are acceler- 
ated by underlying stellar winds or highly coUimated jets, 
the actual mechanism of acceleration remains uncertain. 
The most favoured model, at least for low and interme- 
diate mass stars, is that the molecules are accelerated 
at jet-driven bow shocks. Here we investigate, through 
high resolution numerical simulations, the efficiency of this 
mechanism in accelerating ambient molecular gas without 
causing dissociation. The efficiency of the mechanism is 
found to be surprisingly low suggesting that more momen- 
tunr may be present in the underlying jet than previously 
thought. We also compare the monrentum transferring ef- 
ficiencies of pulsed versus steady jets. We find that pulsed 
jets, and the corresponding steady jet with the same aver- 
age velocity, transfer virtually the same momentum to the 
ambient gas. The additional momentum ejected sideways 
from the jet beam in the case of the pulsed jet only serves 
to accelerate post-shock jet gas which forms a, largely 
atomic, sheath around the jet beam. 

For both the steady and pulsing jets, we find a power 
law relationship between mass and velocity {m{v) ex v~'^) 
which is similar to what is observed. We also find that 
increasing the molecular fraction in the jet decreases 7 as 
one might expect. We reproduce the so-called Hubble law 
for molecular outfiows and show that it is almost certainly 
a local effect in the presence of a bow shock. 

Finally, we present a simple way of overcoming the 
numerical problem of negative pressures while still main- 
taining overall conservation of energy. 

Key words: hydrodynamics - shock waves - ISM:jets and 
outfiows - ISM:molecules 



1. Introduction 

It has been proposed that molecular outfiows, at least from 
low and intermediate mass young stars, may be driven by 
highly coUimated jets (see, for example, Padman, Bence 
& Richer 1997). Although the likely mechanism by which 
such jets transfer their momentum to the ambient medium 
remains unknown, a number of ideas have been put for- 
ward (for a review o f mod els the reader is referred to 
Cabrit, Raga & Gueth 1997 ). Of these the most promising 
seems to be the so-called "prompt entrainment" mecha- 
nism. According to this model the bulk of the molecu- 
lar outfiow is accelerated ambient gas near the head of 
the jet or more precisely along the wings of its associated 
bow shock. Observational support for prompt entrainment 
comes from the spatial coincidence of shocked molecular 
hydrogen bows with peaks in the CO outflow emission 
(e.g. Davis & Eisloffel [l99^ ) . 

Smith, Suttner & Yorke (|1997D and Suttner et al. 



(1997) have carried out a number of 3-D simulations of 
dense molecular jets propagating into a dense medium in 
order to test the prompt entrainment hypothesis. These 
authors found that their simulations reproduced many of 
the observational characteristics of molecular flows includ- 
ing t he so-called 'Hubble law' (see, e.g. Padman et al. 
|l997 ) and stro ng for ward, as opposed to sideways, motion 
(Lada & Fich 1996| ). While such results are encouraging 
for jet-driven models, it is still fair to say that no indi- 
vidual model has yet been able t o plau sibly account for 
all the observations (Lada & Fich 1996 ). Moreover, alter- 
natives to the jet model may be better at explaining the 
observational characteristics of some molecular flows (e.g. 
Padman et al. |1997| and Cabrit et al. |1997|). 



The limited resolution of the 3-D jet simulations of 
Smith et al. ( |1997D , and Suttner et al. I \l997\j , along with 
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the high densities used by these authors, meant that they 
could not resolve the post-shock cooling regions in the 
flow. In addition it was not possible to explore parameter 
space as only a few such simulations could be performed. 
Here we take a somewhat different approach by assum- 
ing low density atomic/molecular jet mixtures, a low den- 
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sity ambient medium and cylindrical symmetry. Although 
such an approach obviously has it limitations, it does allow 
us to explore parameter space more fully and to resolve 
post-shock cooling regions (this might be important, for 
example, if one is to gauge the importance of certain in- 
stabilities). The primary goal of this work is to investigate 
the efficiency of YSO jets in accelerating ambient molec- 
ular gas without causing dissociation of its molecules. 

An additional question we address in this paper is 
whether velocity variations (pulsing) of the jet might en- 
hance transfer of momentum from the jet to its surround- 
ings and thus help to accelerate ambient gas. Pulsing in- 
duces internal shocks which can squeeze jet gas sideways 
(Raga et al. 1993 ). This gas does not interact with the am- 
bient medium directly, but is instead squirted into the co- 
coon of processed (post-shock) jet gas, which separates the 
jet from th e "shr oud" of post-shock ambient gas. Chernin 
& Masson ( 1995 ) however argue that, through the cocoon, 
momentum from the jet may be continuously coupled to 
the ambient flow. 

The properties of the simulated systems in which we 
are interested are as follows: 

— How much momentum is transferred to the ambient 
molecules? 

— Is there a power-law relationship predicted between 
mass in the molecular flow and velocity? 

— What are the proper motions of the molecular 'knots', 
and how does their emission behave with time? 

— Is the so-called 'Hubble law' of molecular outflows re- 
produced under reasonable conditions? 

— Is there extra entrainment of ambient gas along the jet 
due to velocity variations? 

We will discuss each of these points in turn when pre- 
senting our results. 

Our numerical model is presented in §2 and our re- 
sults in §3. Conclusions from this work are presented in 
§4 and a simple way of overcoming the numerical problem 
of negative pressures while still maintaining overall energy 
conservation is given in the Appendix. 

2. Numerical model 

2.1. Equations and numerical method 
The equations solved are 



dt 



= -V • [riYiu) + 2nu2nHk{T) 



|f^-v,.„, 



dt 



-V • [puu + PI] 



1^ = -V-[{e + P)u]-L 



dnux 

dt 
dnu.2 

dt 



= -V • (nn^w) - nH2riHfc(T) 



(1) 
(2) 
(3) 
(4) 
(5) 



^ = - V . ipru) 



(6) 
(7) 



where p, u, P, e and / are the mass density, velocity, 
pressure, total energy density and identity matrix respec- 
tively. Tifj and rifja are the number densities of atomic and 
molecular hydrogen, x is the ionization fraction of atomic 
hydrogen, T is the temperature, J{x,nii,T) is the ioniza- 
tion/recombination rate of atomic hydrogen, k{T) is the 
dissociation coefficient of molecular hydrogen, and r is a 
passive scalar which is used to track the jet gas. We also 
have the definitions 



1 , ^^ u 
e = —pu ■ u H P 

2 Kb 

L = Lrad + EiJ{x, Uh, T) + Ej,k{T) 



(8) 
(9) 



where c^ is the specific heat at constant volume, fee is 
Boltzmann's constant, Ei is the ionization energy of hy- 
drogen and Ed is the dissociation energy of H2. So L is 
a function which denotes the energy loss and gain due to 
radiative and chemical processes, irad is the loss due to 
radiative transitions and is made up of a functio n for losses 
due to atomic transitions (Sutherland & Dopita |l993 ), and 
one f or losses due to molecular transitions (Lepp & Shall 
1983). The second term in L is the energy dumped into 
ionization of H, and the third is that dumped into dissoci- 
ation of H2. The dissociation coefficient k{T) is obtained 
from Dove & Mandy ( 1986| ) and the ionization rate, J, is 
that used by Falle & Raga ( |1995D . 

These equations are solved in a 2D cylindrically sym- 
metric geometry using a temporally and spatially sec- 
ond order accurate MUSCL scheme (van Leer 1977; Falle 



1991). The code uses a linear Riemann solver except where 



the resolved pressure differs from either the left or right 
state at the cell interface by greater than 10% where it uses 
a non- linear solver (following Falle 1996, private commu- 
nication). Non- linear Riemann solvers allow correct treat- 
ment of shocks and rarefactions without artificial viscos- 
ity or entropy fixes. Applying them only in non-smooth 
regions of the flow means that, while the benefits are the 
same, the computational overhead is minimised. This code 
i s an u pdated version of that described in Downes & Ray 
( |1998|) . 

Sometimes negative pressures are predicted by simula- 
tions involving radiative cooling. Typically these are over- 
come by simply resetting the calculated pressure to an 
arbitrary, but small, positive value. However, this involves 
injecting internal energy into the system and this is unde- 
sirable. A fairly reliable way of overcoming this problem 
is discussed in Appendix A. 



2.2. Initial conditions 

Initially the ambient density and pressure on the grid are 
uniform and defined so that the ambient temperature on 
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the grid is 10^ K. The jet temperature is set to lO'^ K. The 
function L is set to zero below this latter temperature as 
the data used in the cooling functions becomes unreliable 
and cooling below this temperature is not dynamically 
significant anyway. In most cases the ratio of jet density 
to ambient density (= 77) is set to 1 (see Table m. The 
ratio — ^ = 9 both inside and outside the jet, unless oth- 
erwise indicated (again, see Table |l|). In all cases the gas 
is assumed to be one of solar abundances. The boundary 
conditions are reflecting on r = (i.e. the jet axis) and on 
z = except where the jet enters, and gradient zero on ev- 
ery other boundary. The computational domain measures 
1500 X 300 cells (but larger in the 77 = 10 simulations), 
with a spacing of 1 x lO^'' cm. We find that the efficiency 
of momentum transfer is sensitive to the grid spacing. We 
performed a number of simulations with different spacings 
and concluded that this is the absolute minimum neces- 
sary to get reliable results. This length should be reduced 
with increasing density. Incidentally this means that ex- 
amining this property at the densities used by, for exam- 
ple. Smith et al. (1997) is impractical. 

The jet enters the grid at = and r < R and the 
boundary conditions are set to force inflow with the jet 
parameters. R is set at 5 x 10^^ cm or 50 grid cells. The 
jet velocity is given by 



4 



(10) 



with Wo ~ 215 km s~^ corresponding to a Mach number 
of 65 and loj are chosen so that the corresponding periods 
are 5, 10, 20 and 50 yrs. Here vi is effectively an ampli- 
tude for the velocity variations where present. The jet is 
initially given a small shear layer of about 5 cells (5 x lO^'' 
cm) in order to avoid numerical problems at the bound- 
ary between the jet and ambient medium. In this layer the 
velocity decays linearly to zero. 

Nine simulations were run with varying values of den- 
sity and velocity perturbation. These are listed in Table 
|l| which also gives the key we will use to refer to the sim- 
ulations. In addition simulations of a purely atomic jet, 
and of a jet with a wide shear layer of almost 25 cells 
(2.5 X 10^^ cm) were run. The velocity of the jet with the 
wide shear layer is given by 



Viet{t,r) 



Wjot(t) 



1 — tanh 



R 



IQi^cm 



(11) 



where v^etit) is given by Eq. |l^. Each system was simu- 
lated to an age of 300 yrs. Although this is very young 
compared to the observed age of stellar jets, it was felt 
that the qualitative behaviour of the system at longer 
times could reliably be inferred from these results. The 
densities chosen are rather low to ensure adequate res- 
olution of the system as described above. Unfortunately 
this precludes the use of the data of McKee et al. ( 1982|) 



Label 


Jet density (cm ^) 


( — ) 


Notes 


A 


10 





- 


B 


10 


0.6 


- 


C 


100 





- 


D 


100 


0.1 


- 


E 


100 


0.2 


- 


F 


100 


0.4 


- 


G 


100 


0.6 


- 


Gl 


100 


0.6 


Wide shear layer 


G2 


100 


0.6 


Atomic jet 


H 


100 


0.0 


77= 10 


I 


100 


0.6 


77= 10 



Table 1. A list of the simulations performed in this work. 
Unless otherwise stated ?7 = 1 



for calculations of the emissions from CO as their calcu- 
lations are only valid for gases of much higher densities. 
As a result we present our findings in terms of the mass 
of molecular gas rather than its luminosity. 

3. Results 

Fig. |l| shows plots of the distribution of number density 
for simulations C and G. The cocoon of the varying jet has 
many bow-shaped shocks travelling through it as a result 
of the internal working surfaces in the jet forcing gas and 
momentum out of the jet beam. It is interesting to note 
that the bow shock of the steady jet is more irregular 
than that of the pulsed jet. This irregularity is probably 
due to the Vishniac instability (e.g. Dgani et al. 1996| ) 
growing at the head of the jet. Presumably the variation 
of the conditions at the head of the varying jet damp- 
ens the growth of this instability. It should be noted that 
the enforced axial symmetry in these calculations makes 
the bow shock appear more smooth than it would in 3D 
simulations. However, the phenomenon of the bow shock 
breaking up occurs in both 2D and 3D simulations. 

We will discuss each of the properties mentioned in §|| 
below. 



3.1. Momentum transfer 

We make use of the jet tracer r to track how much mo- 
mentum has been transferred from the jet to the ambient 
medium. The fraction of momentum transferred from jet 
gas to ambient molecules is 



•^H, 



Z) I , j ( 1 ^ "^y- ) I ^ii I '^H2 "H2 , y- dVij 
J2ij \utj\ntj <m> dVij 



(12) 



where TT-Ha.ij a-nd riy are the number density of molecular 
hydrogen and the total number density in cell ij respec- 
tively, i and j are cell indices in the z and R directions, 
and dVij is the volume of cell ij. This equation is valid 
since the only momentum on the grid originated in the jet, 
and since the simulations are stopped before any gas flows 
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Fig. 1. Log-scale plots of the distribution of number density for (a), simulation C, and (b), simulation G at t = 300 
yrs. The scales are in units of cm^^ 



off the grid. Note that we only consider momentum trans- 
ferred to ambient molecules because we are only interested 
in how efhcient YSO jets are at accelerating molecules, not 
atoms. Thus we ignore ambient molecules which have been 
dissociated in the acceleration process. Table g shows the 
fraction of momentum in ambient molecules for all the 
simulations after 300 yrs. For completeness we also show 
•^totab the total fraction of momentum transferred to am- 
bient gas, whether molecular or atomic. 

The first interesting point to note from Table || is that 
the amount of momentum residing in ambient molecules in 
these simulations is typically an order of magnitude less 
than the total momentum contained on the grid. Note, 
however, how significant amounts of momentum are trans- 
ferred to the ambient medium as a whole, especially in 
those cases where the jet density matches that of its en- 
vironment. This is as one would expect. What is perhaps 
surprising at first is the low efficiency of momentum trans- 
fer to ambient gas that remains in molecular form in the 
post-bow shock zone. 

Comparison between .Ftotai and .Fua for models A and 
C and models B and G clearly shows that the momentum 
transfer efficiency from the jet to the ambient medium 
decreases with increasing density. This result is particu- 
larly marked in the case of post-shock ambient molecules. 
Although further simulations should be performed to con- 
firm this finding, it is physically plausible. Cooling causes 
the bow shock to be narrower (i.e. more aerodynamic) 



than in the adiabatic case, thus reducing its cross-sectional 
area. Obviously this leads to a reduction in rate at which 
momentum is transferred from the jet to its surroundings. 
The fact that the effect is more marked for post-shock 
ambient molecules must reflect changes in the shape of 
the bow (as opposed to pure changes in its cross sectional 
area) with increased cooling. 

Since we are simulating systems here which are prob- 
ably of low density in comparison to typical YSO jets, 
our results suggest that radiative bow shocks, from at 
least heavy and equal density jets (with respect to the 
environment), are not very good at accelerating ambient 
molecules without causing dissociation. This result also 
points to the fact that in the case of such jets, the jet may 
carry much more momentum than one might naively es- 
timate based on a rough balance with the momentum in 
any associated observed molecular flow. 

We now turn to differences in the efficiency of mo- 
mentum transfer in pulsed versus steady jets. The topic 
of differences in entrainment rates will be discussed more 
fully in §3.5. Fig. shows grey-scale plots of the distribu- 
tion of |m| and of jet gas for models C (^ = 0) and G 
(it ^ 0.6). Comparison between the steady jet and the 
varying velocity jet suggests that momentum is indeed 
being forced out of the beam of the varying jet by the in- 
ternal worki ng surfaces as predicted for example by Raga 
et al. ( 1993 ). It is also interesting to note the similarity 
between the distribution of velocity and the distribution 
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Model 


■^Ha 


•>^total 


7 


A 


0.21 


0.58 


2.42 


B 


0.20 


0.55 


2.93 


C 


0.08 


0.36 


2.37 


D 


0.07 


0.36 


2.08 


E 


0.07 


0.36 


1.81 


F 


0.09 


0.33 


2.31 


G 


0.10 


0.38 


2.98 


Gl 


0.10 


0.39 


2.02 


G2 


0.10 


0.39 


3.75 


H 


0.06 


0.16 


1.58 


I 


0.08 


0.21 


2.44 



Table 2. The proportion of momentum on the grid re- 
siding in ambient molecules (^Hs) and in all ambient gas 
(-^totai) at i = 300 yrs. Also given is the value of 7, the 
power-law index for the mass- velocity relationship 



of jet gas in both simulations. Moreover it is clear that the 
momentum leaving the jet beam is dumped into jet gas 
which has been processed through the jet-shock and inter- 
nal working surfaces and now forms a cocoon around the 
jet itself. Since this gas is largely atomic (most of it having 
passed through strong shocks), this effect does not directly 
lead to extra acceleration of molecular gas. However, the 
ejected momentum could conceiveably pass through the 
cocoon of jet gas eventually and go on to accelerate am- 
bient molecules. The wings of the shocks caused by the 
internal working surfaces in the jet have encountered the 
edge of the cocoon by the end of these simulations but, 
even so, the fraction of momentum in ambient molecules 
varies by not more than 3% as a result of th e velo city vari- 
ations. This result was noted by Downes ( |1996 ) for slab 
symmetric jets, but here we extend this result to cylindri- 
cal jets with a variety of strengths of velocity variations. It 
is also interesting to note that, from comparisons between 
simulations G and I, the efficiency of momentum transfer 
(in particular to molecular gas) is not very sensitive to ry. 



3.2. The mass-velocity relationship 

We do find a power-law relationship between mass of 
molecular gas and velocity. If we write 



m{v) oc f ''' 



(13) 



then we find that 7 lies between 1.58 and 3.75 and that 
7 tends to increase with time, in agreement with Smith 
et al. ( |1997D . Table | shows the values of 7 at t == 300 
yrs for all the models. These values are consistent with 
observations (e.g. Davis et al. 1998), and also with the 



analytical model presented in the appendix of Smith et al. 
( |1997D for the variations of mass with velocity. However, it 
is important to emphasise that what is actually observed 
is a variation in CO line intensity with velocity. CO line 
intensity is directly proportional to mass, in the relevant 
velocity channel, providing we are in the optically thin 



Velocity (km/s) 

Fig. 3. Plot of the relationship between the molecular 
mass and velocity for simulation G at t = 300 yrs as- 
suming the jet moves at an angle of 60° to the plane of 
the sky. Note how a power-law (dashed line) fits the data 
quite well 



regime and the temperature of the gas is higher than the 
excitation temperature of the line (see, e.g. McKee et al. 



19821 ). Note that Smith et al. ( |l997l ) incorrectly state that 
the channel line brightness scales with v^dm{v). Fig. ^ 
shows a sample plot of the molecular mass versus velocity 
for the jet moving at an angle of 60° to the plane of the 
sky. We do not see the jet contribution in the velocity 
range chosen here. 

The molecular fraction in the jet has a marked influ- 
ence on the value of 7 predicted by these models as we can 
see by comparing the results for simulations G and G2. In 
fact, 7 increases with decreasing molecular abundance in 
the jet. This is due to the reduction in strength of the high 
velocity jet component. 

It appears that 7 does not depend in a systematic way 
on the amplitude of the velocity variations. Note also that 
the introduction of a wide shear layer dramatically reduces 
7. This is due to the fact that more gas is ejected out of the 
jet beam (because of the more strongly paraboloid shape 
of the internal working surfaces) and this accelerates the 
cocoon gas, leading to a stronger high velocity component. 
In addition, a wide shear layer causes the bow shock to 
be more blunt. It can be seen from the analytic model of 
Smith et al. ( |1997D that this also leads to a lower value of 
7. It is interesting to speculate that lower values of gamma, 
which may be more common in molecula r outfl ows from 
lower luminosity sources (see Davis et al. 1998 ) could re- 
sult from such flows having a higher molecular fraction in 
their jets and perhaps a wide shear layer. 

The behaviour of 7 with viewing angle is the same as 



that noted by Smith et al. (1997). The actual values of 7 



obtained by these authors are somewhat lower than those 
obtained here. However, since our initial conditions are so 
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Fig. 2. Plots (a) and (b) show the distribution of the jet tracer variable and |m| for simulation C respectively. Plots 
(c) and (d) show the same quantities for simulation G. Plots (a) and (c) use a linear scale and a value of 1 indicates 
pure jet gas. Plots (b) and (d) also use a linear scale and are in units of km s~^ 



different, and since 7 is dependent on the shape of the bow 
shock, this discrepancy is not disturbing. 



3.3. H2 proper motions and emissions 

We measured the apparent motion of the emission from 
the internal working surfaces. Near the axis of the jet this 
emission moves with the average jet speed (i.e. vq), as 
would be expected from momentum balance arguments. 
However, there are knots of emission arising from the bow 
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shock itself and these move much more slowly (~5-15% of 
the average jet speed) with the faster moving knots being 
closer to the apex of the bow. This is in agreement with 
the observations of Micono et al. ( 1998| ). 

Fig. ffl shows the emission from the S(l)l-0 line of H2 
for model G. There is little emission from the cocoon since 
the cocoon gas has been strongly shocked in the jet shock 
and so is mostly atomic. We can also see that the emission 
becomes more intense as we move away from the apex of 
the b ow sh ock, as reported by many authors (e.g. Eisloffel 
et al. 1994). It is also clear that the internal working sur- 



faces in the jet are giving rise to emission in this line. We 
can see that the emission begins to die away as we move 
away from the jet source. This is in agreement with obser- 
vations of, for example, HH 46/47 (Eisloffel et al. [1994D 
where the emission from the knots appears close to the jet 
source and then fades away. 

This decrease in emission happens for two reasons. 
The first is that the shocks in the jet become weaker as 
they move away from the source simply because the ve- 
locity variations, which give rise to the shocks in the first 
place, are smoothed out by the shocks (see, for example 
Whitham (1974). In addition, the mass flux through an 



individual shock decreases with time due to the divergent 
nature of the flow ahead of each internal working surface. 
This means that the emission will decrease because there 
is less gas being heated by the shock. 

3.4. The 'Hubble law' 

We have found that the so-called 'Hubble law' (e.g. Lada 
& Fich 1996) is reproduced in these simulations. Fig. S 
shows a position velocity diagram calculated from simu- 
lation G assuming that the jet makes an angle of 60° to 
the plane of the sky. This diagram is based on the mass of 
H2 rather than intensity of CO emission. There is a grad- 
ual, virtually monotonic, rise in the maximum velocity. It 
is also worth noting that near the apex of the bow shock 
the rise in the maximum velocity present becomes steeper. 
These properties are related to the shape of the bow shock 
as gas near the apex of the shock is moving away from the 
jet axis at higher speed than that far from the apex. 

As a very basic model of this, suppose we represent 
the contact discontinuity between the post-shock jet and 
ambient gas to be an impermeable body moving with ve- 
locity V through a fluid whose streamlines will follow the 
surface of the body. See Fig. o for a schematic diagram of 
the system. Let this surface be described by the equation 



z = a — r 



(14) 



where s > 2 and a is the position of the apex of the bow 
shock on the z axis (see, e.g.. Smith et al. 1997 ). Since the 
contact discontinuity is a streamline of the flow we get 



that the ratio of the z-component to the r-component of 
the velocity is simply 



n 



= s la — z\ 



(15) 



Note that this is the negative of the slope of the bow shock. 
This is because of our choice of the bow shock pointing 
to the right, and hence the z component of the velocity 
will be negative. If we assume the post-shock velocity to 
be vi{z) (related to v by the shock jump conditions), it is 
simple to show that 



Vriz) = 



Vi 



ViTW 



(16) 



Finally, after some simple algebra, we can write down the 
velocity along the line of sight as a function of z by 



V\os{z) 



Vl 



VTT^ 



{cos a + 7?. sin a} 



(17) 



where a is the angle the bow shock makes to the plane of 
the sky. In fact we can write down vi{z) if we assume that 
the bow shock is a strong shock everywhere, thus yielding 
a compression ratio of 4 (from the shock conditions) . It is 
easy to derive that 



vi{z) 



V cos (arctan(7?.)) \j l-TZ"^ 

" 16 



(18) 



This formula yields a shape for the position-velocity 
diagram which suggests the Hubble-law and is similar to 
diagrams generated from these simulations assuming that 
the flow is in the plane of the sky. This indicates that the 
'Hubble-law' effect is, at least partly, an artifact of the 
geometry of the bow shock. 

3.5. Entrainment 

As is clear from Fig. g there is not much extra entrainment 
of ambient gas resulting from the velocity variations. If the 
velocity variations were to involve the jet 'switching off' 
for a time comparable to the sound crossing time of the 
cocoon then we would expect ambient gas to move toward 
the jet axis and probably be driven into the cocoon when 
the jet switches on again. This does not happen in these 
simulations where the maximum period of the variations 
is 50 yrs and the amplitude is at most 60% of the jet 
velocity. However, there is a small amount of acceleration 
of molecular gas at the left-hand boundary of the grid. 
This effect is quite small, but may grow over time. It is 
also possible, however, that this effect is due simply to the 
reflecting boundary conditions. 

These simulations show very little mixing between jet 
and ambient gas except very close to the apex of the bow 
shock. This means that the high velocity component of CO 
outflo ws often observed along the main lobe axis (Bachiller 
1996 ) is difficult to explain without invoking the presence 
of CO gas in the jet beam itself just after coUimation. 
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Fig. 4. Log-scale plot of the distribution of emission from the S(l)l-0 2.12/i line of H2 for simulation G after 300 yrs. 
The scale is in units of erg cm""^ s~^ 
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Fig. 5. Contour plot of the position- velocity diagram (of 
mass rather than emission) for simulation G after 300 yrs 
assuming the jet moves at an angle of 60° to the plane of 
the sky. Note the gradual rise of the maximum velocity as 
we move away from the source. The contours are logarith- 
mic running from 3 x 10^^ to 3 x 10^^ g. The contribution 
from the jet is removed from this diagram to make the 
effect clearer 



4. Conclusions 

It is now generally agreed that low velocity molecular out- 
flows represent ambient gas that is somehow accelerated 
by highly coUimated and partly ionized jets, at least in 
the case of low and intermediate mass YSOs. Most of this 
acceleration is thought to be achieved at the head of the 
jet through the leading bow shock (the so-called prompt 
entrainment mechanism). In this paper we have examined 
through many axially symmetric simulations the efficiency 
of the prompt entrainment mechanism as a means of trans- 
ferring momentum to ambient molecular gas without caus- 
ing dissociation. It is found, as one would expect, that the 
fraction of jet momentum transferred to the ambient envi- 




Fig. 6. Schematic diagram of the setup used to derive the 
Hubble law for position- velocity diagrams. The fluid mo- 
tion is shown by vectors. After contact with the bow shock, 
the fluid is assumed to flow along the contact discontinuity 
with a smaller, but constant, velocity. See text 



ronment depends on the jet/ambient density ratio. More 
importantly, we see that cooling, which is particularly im- 
portant at higher densities, decreases the fractional jet mo- 
mentum that goes into ambient molecules. It would seem 
on the basis of the simulations presented here that both 
heavy and equal density (with respect to the environment) 
jets with radiative cooling have very low efficiencies at ac- 
celerating ambient molecules without causing dissociation. 
In part this is because cooled jets have more aerodynamic 
bow shocks than the corresponding adiabatic ones (i.e. 
they present a smaller cross sectional area to the ambient 
medium). The actual shape, however, of the bow shock 
also seems to to be important as the decrease in momen- 
tum transferred to the ambient medium seems to affect 
the acceleration of molecules more so than atoms/ions. 
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We have also tested whether pulsed jets are raore ef- 
ficient at transferring momentum to the ambient medium 
than the corresponding steady jet with the same aver- 
age velocity. Somewhat surprisingly we found that even 
relatively large velocity variations do not give rise to sig- 
nificant changes in the amount of momentum being de- 
posited in ambient gas. Fundamentally this is because in 
high Mach number jets, even with cooling, there is very 
little couphng between the jet's cocoon and the "sheath" 
(i.e. the post-shock ambient gas). This lack of coupling is 
also the reason why turbulent entrainment is not signifi- 
cant in YSO jets. 

Our simulations were also used to model the expected 
variation of mass with velocity in molecular flows. We 
found relatively large 7 values, i.e. mass should decline 
steeply with velocity, in line with observations. Interest- 
ingly we found that a wide shear layer and an increas- 
ing molecular component in the jet reduced 7. If, as one 
might expect, such conditions are common among out- 
flows from low luminosity YSOs, this could explain their 
observed lower values for 7. Finally we have shown that 
the so-called Hubble Law for molecular outflows is almost 
certainly a local effect in the vicinity of a bow shock. 

Appendix A: Overcoming negative pressures 



1 ?i+ — 
fluxes F . , 1^ instead of the second order ones. Note that 

it is necessary to adjust the neighbouring cells also so that 
overall conservation is maintained. 

While this fix does not work for all problems, it was 
found to be rather effective in the simulations presented 
here where very strong rarefactions are produced both at 
the edge of the jet at the z — boundary, and within 
the jet itself where the enforced velocity variations can 
cause problems. This fix can be implemented very easily 
by ensuring that the flux vectors around problem cells are 
stored. It has the advantage that it does not involve losing 
conservation of any of the physically conserved quantities. 
Reducing the scheme to first order in certain regions of 
the grid is not a significant problem because this happens 
around shocks anyway in order to maintain monotonicity. 
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As noted in §2.1 it is common for conservative numerical Kererences 



codes to predict negative pressures under certain condi- 
tions. Flows giving rise to such problems are usually highly 
supersonic, and diverging. The difficulties occur because 
the ratio of internal to total energy goes like -p^, where 
M is the Mach number of the flow. Therefore, in a highly 
supersonic flow, the fractional error required to predict a 
negative pressure is rather small. The introduction of en- 
ergy losses exacerbates this problem. This is due to the 
fact that relatively weakly diverging flows, for example, 
can become supersonically diverging when the system is 
cooled because the sound speed is reduced. 

Schemes which are second order accurate in space tend 
to produce more negative pressures than first order ones. 
This is because first order schemes dissipate strong fea- 
tures quickly so that strongly diverging flows rarely oc- 
cur. It seems reasonable, then, to invoke a scheme which 
is first order in space whenever a negative pressure is pro- 
duced as this will introduce extra dissipation. This extra 
dissipation may eliminate the negative pressure by moving 
some extra energy from neighbouring cells into the prob- 
lem one. The scheme used by the code in this work can be 
summarised as follows: 



jjn+l = [/" A 



-,n+T 



i+i 



(A.l) 



where C/" and '^Ff are the state vector and second order 
flux calculated at t = nAt and x = iAx respectively, and 
A = ■^. If this scheme produce s a negative pressure at 



1 then simply apply Eq. |A.l| using the first order 
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